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INTRODUCTION 


Under the present Air Force Contract F19628-81-K-0012, 
'Ground Response In Alluvial Basins Due to Seismic 
Disturbances", the accomplishments of the contract are: 

(I) Five preliminary basin models for SH-waves and three 
for elastic (P and SV) waves have been studied using the 
finite element method. The results of these problems are 
presented either by the conventional synthetic seismograms or 
by tha detailed snapshots of the wave fields. These results 
should give the effects of the basin and range like structure 
on wc.'.'O scattering, diffraction, and propagation. 

(II) The problem of undesired reflections and diffractions 
from artificially terminated boundaries In the finite element 
modeling has been re-1nvestIgated. A non-reflect1ng-boundary 
Is implemented. By using this type of boundary, artificial 
reflections can be eliminated completely, before a given wave 
Is reflected at that same boundary more than once. 

(Ill) A new finite element technique has been successfully 
developed to simulate different types of earthquake and 
explosive source mechanisms, including the following cases: 

(1) One Degree of Freedom SH-wave Problems with: 








2 


(a) a concentrated line source, 

(b) a concentrated coupled line source. 

(2) Two Degrees of Freedom Elastic-Wave <P and SV) 
Problems with: 

(a) a d1rectIona 1-force line source, 

<b) an omni-directional line source, 

(c) a s1ngle-couple-wlthout-moment line source, 

(d) a sIngle-couple-wlth-moment line source, 

(e) a douple-couple line source. 

(IV) The two-dimensional element-oriented finite element 
codes for the cases of elastic (P and SV-wave) and SH-wave 
have been modified and Improved by: 

(1) Implementing the Effective Excitation (EE) method, 

(2) Introducing the relative coordinates of the nodal 
points, 

(3) Introducing the restart and back-up capability, and 

(4) replacing the 4CST(Four Constant Strain Triangles) 
algorithm Instead of using the Averaging 2CST 

a I gor 1thm. 

(V) The two-dimensional noda1-po1nt-or1ented finite 








element codes for the case of elastic waves has been 


developed. By using this new version, a large amount of 
computing time can be saved. 

(VI) The element-oriented finite element codes for the 
cases of elastic waves and SH wave, which were originally 
written for PRIME 750, have been adapted to VAX 11/780 at 
AFGL . 
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SCIENTIFIC ACCOMPLISHMENTS 


I. FINITE ELEMENT MODELING OF BASINS WITH DISPLAY OF 


SYNTHETIC SEISMOGRAMS 


(1) SIMPLE MODELS FOR SH WAVES 

The following three simple basin models for SH waves 
have been investigated. 

(a) Model l: 

At the suggestion of Mr. J. Battls at AFGL, we first 
model an homogeneous truncated medium as shown In Figure 1, 
with 500 quadrilateral finite elements. The model assumes 
shear velocity v g - 2,900 m/sec, density ■ 2.67 gm/cnr j 
grid size Ax ■ Ay ■ 85m. A uniform SH wave forcing 
function of the first derivative of a Gaussian function 
(Figure 2) Is loaded on the lower boundary AA' of the 
med 1 urn. 

The synthetic seismograms of the displacement u along 
the vertical planes OA and BC, with a sampling space 
x » 425m and time step t ■ 0.02 sec, are shown In 
Figures 3 and 4. These two se1smogram-vert 1ca1 profiles are 
exactly Identical. Figures 5 and 6 show the Identical 
responses along the horizontal surfaces 00' and MM'. As 
expected, all of these responses due to a uniform external 
SH forcing function along AA', are Identical and assure the 
correctness of the finite element results. In Figures 3 and 
4, the arrivals of the direct and reflected pulses from AA 1 


and 00' can be clearly Identified. 


" -V 
















Figure 1. Finite Element Model 1. 
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thetic Seismogram Along Vertical Plane OA 




























igure5. Synthetic Seismogram Along Free Surface 00' For 











igure 6. Synthetic Seismogram Along Horizontal Surface MM' For Model 














(b) Model 2: 


The exact same model as (a) above, except the upper 
left corner c e 25 elements of medium 1 Is replaced by medium 
2, with a shear wave velocity, v g2 * 1,750 m/sec, and 

3 

density p 2 * 2.4 gm/cm (Figure 7). The excitation 

function Is Identical to that In Model 1. 

Figures 8 and 9 show the complex results of the 

displacement uz along the vertical planes OA and BC. The 

Inclusion of medium 2 Introduces a great deal of complexity 

of the waves. In Figures 8 and 9, the arrivals of the 

direct waves and first reflections from the free surface 00' 

are still Identifiable; however, the diffracted waves and 

mu 111-ref 1ected waves In medium 2 give the complex waveforms 

of the reflected waves at 00'; the reflected waves from AA' 

are distorted. As expected, the responses of the 

displacements u„ along any horizontal surfaces are no 

z 

longer Identical as compared with the situation In Figure 6. 
Figures 10 and 11 are the synthetic seismograms along the 
horizontal surfaces 00' and MM'. 














Reflected at 00' /yReflected at 00' 



nthetic Seismogram Along Vertical Plane OA for Model 
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(c ) Model 3: 

In this section, two partial bas-lns are modeled. Both 
basins use the same homogeneous material as Model 1 above, 
but with a different geometry, as shown In Figures 12 
(Model 3-a) and 13 (Model 3-b>. In Figure 12, an SH wave 
line source Is excited at Point S on the free surface of the 
deeper basin, while the receivers are located on the free 
surface of the shallower basin. In the region BC. For Model 
3-b (In Figure 13), an SH source Is located In the shallower 
basin, while the receivers are In the deeper basin. Again, 
the forcing function Is the first derivative of a Gaussian 
function with Its center frequency f c ~ 2.5 Hz. The 
calculated displacements u z for Models 3-a and 3-b 
obtained along ABC In both oodels are shown In Figures 14 
and 15, respect 1vely. In Figures 14 and 15, at the 
corner B, as expected, waves change their magnitude and 
phase. In the case of Model 3-b, the reflection from plane 
DB causes drastic changes In phase and amplitude. From the 
distortion of these waves, shown In Figure 14, we may 
conclude that, since the arrivals at the position ABS’ are 
the diffracted waves generated at 0, the effect of the 
diffraction due to the corner Is small In comparison with 
the arrivals of the direct waves received at the position 
ABS' In Model 3b (Figure 13). 

















Figure 13 . Preliminary Basin Model 3 b 












































































(2) SIMPLE MODELS FOR ELASTIC WAVES (P & SV WAVES) 


(a) Prlllmlnary Basins Model 3-A and Model 3-B 

We have studied two preliminary basin models as shown 
In Figures 16 and 17 for the case of elastic (P and S) 
waves. In Figure 16, the homogeneous medium of the 

preliminary basin model <v = 5 ,01313 m/sec, v * 2,900 m/sec, 

3 P S 

= 2.67 gm/cm ) Is excited by a vertical forcing function 

F Q < t) as shown In Figure 18 at Point S on the free surface 
of the deeper basin DC. The receivers are located on 
the free surface of the shallower basin In the region BC. 
In Figure 18 the source of a vertical forcing funtlon Is 
located In the shallower basin. The receivers are located 
In the deeper basin. Again, the forcing function Is the 
first derivative of a Gaussian function with Its center 
frequency f c ~ 2.7 Hz. Figures 19 and 20 are the results 
of the vertical and horizontal displacement components along 
plane ABC for both Model 3A and Model 3B, showing the 
complexltly of the wave arrivals, respectively. 






















Figure l«. Source Function for'Elastic Cases. 

















Lgure l 9 Synth 

















































( b) Two Ajacent Alluvial Basins for 
Elastic Waves <P and S Waves) 


The finite element model of the two adjacent alluvial 
basins, A and B, Is shown In Figure 21. The basin Is filled 


with 3 km of 


sed1ments 


( v _ a 3,500 m/sec , 

P2 



v * 1,905 m/sec, * 2.5 gm/cm^), which are overlain by a 
surface low velocity layer (0.5 km thick, v ^ ** 2,200 m/sec, 

3 p 

v g ^ * 1,180 m/sec, ■ 2.3 gm/cm ). The bedrock of the 
basement Is assumed to have v ^ = 5,500 m/sec, 

v _ * 3,175 m/sec, o = 2.7 gm/cm^. The basin B Is 15 km In 

s 3 / 3 

width and filled with the same materials as basin A. The 
basement Is bounded by a rigid halfspace. The total length 
of the two alluvial basin model Is 50 km. The thickness of 
the model Is 7 km. A vertical force load Is located In a 
smaller basin (Basin A of Figure 21) with a width of 10 km. 
The source forcing function, F Q , used In the finite element 
calculation Is that of the first time derivative of a 
Gaussian function as shown In Figure 22. In order to obtain 
a stable and accurate solution, we choose a time step, At, 
that satisfies the Courant-F rledr fckson-lewy (CFL) condition 
and an element dimension s to be small enough so that one 
wavelength can be approximately simulated by 10 A s. In 
doing so, the width of the forcing function N At should be 

greater than a factor of about (8.3 v At)/v , where v 

max min max 

and v are the largest and the smallest velocities in the 

min 

whole finite element model. 

For tha finite element alluvial basin model, we use a 
source function having a pulse width of 160 At (N=4J). A 
total of approximately 70,000 elements is required for such 
a model. Using a center frequency of the forcing function 
of about 1 Hz, the total number of equations of motion is 









140,000. Figures 23-A to 24-B are the preliminary results 
obtained from the model. Figures 23-A and 23-B are the 
synthetic seismograms of the vertical and horizontal 
displacements at various nodal points along the surface SO 
(see Figure 21). These d1sp1acaments are normalized to the 
displacement of the forcing function at the loading 


position 


Figure 24-A and 24-B are these synthetic 


seismograms along nodal points BO with a scale 10 times 
larger than that used In Figures 23-A and 23-B. At the 
corners B and C, the arrivals are drastically changed In 
magnitude and phase due to the complexity of the basin 
structure. 














Basins Due to Seismic Disturbances. 
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II. FINITE ELEMENT MODELING OF BASINS WITH DISPLAY OF 
SNAPSHOTS 

In order to obtain a clear dtsplay of the results of 
basin modeling, an Aldridge Laboratory computer code for 
snapshots has been Implemented to our present finite element 
computer programs. These snapshots provide a means to 
visualize the time history of the wave interactions and to 
Identify the resulting wave shapes, the direction of 
propagation, and the nature of mode conversion. For 
demonstrat1 on, snapshots have been made for the following 
three cases: 

Case It Model 3-A for Elastic Waves 

Figures 25-A to 25-G are snapshots of the synthetic 
wavefleld at 0.3 sec Intervals. Vertical and horizontal 
displacement components are given separately. Each snapshot 
Is Independently scaled and the amplitude Is magnified 10 
times. Black Is positive displacement and white Is negative 
displacement. Diffracted waves from the square corners are 
clearly Identified after 0.9 sec. 



























Pattern at 













Case 2: Model 3-A for SH-Wave 


In Figure 26, a SH-wave source , F^(t) , Is located at 
the corner S of the lower basin. The source function Is 
shown In Figure 2. Figures 27-A to 27-0 are the snapshots 
of the synthetic wavefield at 0.2 sec Intervals. In this 
simple one-degree-of-freedom problem, reflections from the 
free surfaces and diffractions from the corners can be 
clearly Identified. In Figures 27-M to 27-0, we observe 
that the waveflelds are contaminated by the undesired 
reflections from the side and the bottom boundaries. 












Vit-w of Wave Pattern at Cl./ sec after 
















v 



the betonrtt ion of the Sotu. ♦- the Detonal ion of the Sourc 











the Detonation of the Source 






















the Detonation of the Soui< 















50 


Case 3: Model 4 for SH-Wave 


Figure 28 Is the geometry of a preliminary basin model 
with an Inclined slope of 45° . A SH-wave source function, 
F Q (t), of the type shown In Figure 18 Is located at the 
corner S . Figures 29-A to 29-0 are the snapshots of the 
synthetic wavefleld at 0.2 sec Intervals. Two seconds after 
detonation of the source undesired reflections appear (after 
Figure 2 9 — L ) as a result of reflections from the left side 
of the artificially terminated boundary and 2.6 seconds 
(after Figure 29-M) due to the bottom boundary. Comparing 
the results of this case with those of Case 2, different 
phenomena can be seen In Figures 29-K to 29-0. Because of 
the sloped structure In the neighborhood of the source 
location, the medium Is quiescent until 2.2 sec after 
detonation of the source. However, In Case 2, as the basin 
wall boundaries Is vertical, there Is no such quiescent 
period. Reflections from the vertical portion of the basin 
free surface wall boundaries appear to cause diffractions as 
soon as the waves arrive at the walls. 















Snapshot View of Have Pattern at 0.4 aec after 







Snapshot View of Wave Pattern at 0.6 aec after I 1 u V c _ I) . a Snapshot View of Wave Pattern at 0.8 aec after 













Snapshot View of Wave Pattern at 1.2 *ec after 
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III. MODIFICATION AND IMPROVEMENT OF THE TWO DIMENSIONAL 

ELEMENT COMPUTER CODES 

(I) NON-REFLECTING BOUNDARIES 

Truncation of finite element models generates 
undesirable reflections from the truncated boundaries. 
Numerical analysis techniques that attempt to eliminate 
these undesirable reflections generally Inlolve either 
viscous boundary or superposition boundaries to simulate the 
transparent boundaries. In frequency-domain finite element 
computer modeling, Lysmer and Kuhlemeyer (1969) and White et 
al. (1977) proposed viscous dashpots to damp out most of 
the reflections. However, within the framework of the time 
domain solutions, the viscous dashpot technique Is good only 
for elastic body waves but not for surface waves. Also, the 
damping Is frequency and incident angle dependent. Smith 
( 1973) proposed a superpos111 on method by adding two 
separate solutions, one with a Dlrlchlet and one with a 
Neumann boundary solution, to eliminate the artificial 
reflections. Although Smith's formulation Is Independent of 
both frequency and Incidence angle, It requires 2 n complete 
dynamic solutions If there are n reflecting surfaces. This 
method also falls when a given wave Is reflected at the same 
boundary more than once. 

After considerable effort, we have succeeded In 
refining Smith's superposition method, which Is referred to 
here as modified Smith superposition method for 
non-reflecting boundaries, and It summarized as follows: 

(A) For one degree of freedom problems (acoustic or SH 
waves ) 

(a) Calculate the solution for Neumann's problem for 
each nodal point of the whole structure. 





m 















(b > Divide the solution obtained from Neumann’s problem 
by 2 at the boundaries to obtain the expected solutions at 
the boundaries of each time step. 

(c) Use these expected values at boundaries as feedback 
boundary conditions to calculate the expected responses of 
the Interior nodal points for the same time step. 

(B) For two degree of freedom problems (two-dimensional 
elastic waves ) 

(a) Calculate the solution for the mixed boundary 
conditions 1, l.e., normal displacement and tangential 
stress at the boundary are zero. 

(b) Calculate the solution for the mixed boundary 
condition 2. l.e., tangential displacement and normal stress 
at the boundary are zero. 

(c) Divide the tangential displacement at the boundary 
obtained from (1), and the normal displacement at the 
boundary obtained from (11) by 2 to obtain the recovered 
displacement components at the boundaries at each time step. 

(d) Use these expected values of the displacements at 
the boundaries as the feedback boundary conditions to 
calculate the expected responses of the Interior nodal 
points for the same time step. 

The proposed approach avoids the need for 2 n complete 
solutions if there are n reflecting surfaces. However, this 
method still falls when a given wave Is reflected at the 
same boundary more than once. 

The following are two finite element moda1s to 


demonstrate the cancellation of reflections. 






Example 1 : 

An infinite long elastic plate, bounded by an upper 
free surface and a lower rigid surface, with wave velocities 
Vp ■ 2,000 m/sec. Vs = 1,155 m/sec, density 

* 2.4 gm/cm ^ (Figure 30a > . Symmetry Is used by letting 
the horizontal displacement along z-axts be zero to reduce 
the size of the problem (Figure 30b). We use an Internal 
forcing function, F , with a time step = 5.0x10 ^ sec. The 
observation points are located at A and B. The grid sizes 
are Ax * 16 m and Az ”15 m. First, we construct a model of 
sufficiently large size (150 x 46 elements), to obtain the 
required solution before the arrivals of the reflection from 
the side boundaries (MM' and NN' In Figure 30a>. Then we 
construct a model of smaller size (20 x 46 elements) with 
the application of above described non-ref 1ec11ng boundary 
conditions along BB 1 to simulate an Infinite medium In the 
x-dlrectlon. Figures 31-A and 31-B compare the displacement 
components observed at point B between the two different 
kinds of the solutions. The dotted curves are the expected 
results; the solid curves are the averaged values obtained 
woth free and rigid boundary conditions on BB'. The 
deviation between the two different kinds of solution starts 
at the twenty-first time step numerically. Figures 32-A and 
32-B show the comparison of the expected results and the 
non-ref 1ect1ng boundary solutions based on the modified 
Smith superposition method observed as the same boundary 
point , B. Exact cancellation of the reflected waves at the 
boundary has been effectively achieved until the boundary 
BB' becomes the second-time reflecting boundary for a given 
ray. Figures 33-A and 33-B are the comparison of the two 
kinds of solutions observed at point A; the deviation would 
have started at the thirty-first time step without applying 
the non-reflecting boundary conditions. 









FINITE ELEMENT MODEL 
















FINITE ELEMENT MODEL 
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Figure 31b. Comparison of Horzontal Displacements. 
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Figure 32b. Comparison of Horizontal Displacements 
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OBSERVED AT POINT A (MODEL I ) 



Figure 33b. Comparison of Horizontal Displacements. 
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Example 2: 

In Example 1 we demonstrate the effect of using 
non-ref 1ec11ng boundary conditions on a single boundary. In 
this second example, the non-ref 1ect 1ng boubdary conditions 
are applied to two boundaries. The model Is a half-space 
with the same elastic medium as Example 1 with an Internal 
source (Figure 34a). Again only half of the problem Is 
modeled due to symmetry (Figure 34b). Two models of 
different sizes are constructed: the large one (60x60 
elements) generates the expected solutions; the small one 
(20x20 elements) Is to test the non-ref 1ect1ng effect. The 
small model contains only OABC. Figures 35a and 35b show 
the responses as observed at corner C(see Figure 34). The 
expected responses are obtained by using the large model 
containing OA'B 1 C', and the other solutions are obtained by 
using the small model In which both AC and BC are Imposed 
with free boundary conditions. It Is evident that the 
expe 'ed solutions for both the vertical and horizontal 
components do not agree with the traction-free boundary 
solutions. Analytically, corner C Is singular. However, In 
the finite element algorithm, the non-ref 1ec11ng boundary 
conditions based on the modified Smith superposition method 
can be successfully applied to point C. Figures 36a and 36b 
show the agreement between the above expected solutions and 
the non-ref1ec11ng boundary solutions. Again, this 
non-reflecting technique falls when a given wave Is 
reflected at either boundary more than once. 
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FINITE ELEMENT MODEL II 


Figure 34a. Finite Element Model for an Elastic 

Half-Space. (Region A is the simulated 
half-space with the use of non-reflecting 
boundaries.) 
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Figure 35b. Comparison of Horizontal Displacements. 























































-* TIME STEPS ( = 5 msec) 

Figure 3ba. Comparison of Vertical Displacements. 





























































( 2 > EFFECTIVE EXCITATION (EE) METHOD 


As mentioned In the Final Report of Contract F49620-77 
-C-0130, "Elastic and Visco-elastic Wave Scattering and 
Diffraction" (1977-1980), we used the direct step-by-step 
explicit centra1-d1fference time Integration scheme In 
AFEA3. Consequently, when the structure Is subjected to an 
external disturbance, only the neighboring nodal points of 
the source are primarily excited while all other nodal 
points remain virtually undisturbed. The region of 
excitation motion expands with Increasing time steps. For 
Instance, we first consider a two-dimensional finite element 
mesh as shown In Figure 37. If the source Is located at 
point S, for the first time step we need only M> to 
calculate the displacements, velocities, accelerations of 
those nodal points within the area A^B^C-^D^, Including the 
nodal points, and only (11) to assemble the stiffness and 
mass matrices for those elements within the area ^2^2^2• 
For the second time step, only those nodal points within the 
areas and ^2^2 ^2^2 anc * e ^ ement - s vvlthln the area 
A^B^C-jD^ are effectively excited, and so on. The same 
concept can o» genera 1 (zed to the three-dimensional case. 
In the Aldridge preliminary two-dimensional finite element 
code AFEA, the displacements, velocities, accelerations for 
every nodal point must be calculated for every time step. 
Consequently, the EE-method saves a large amount of 
computing time In the early stage of computation 
particularly for large finite element models. 


















(3) INTRODUCING RELATIVE COORDINATED OF THE NODAL POINTS 


In the old version of the AFEA, the mesh generator 

generates the value of the coordinates of the four nodal 
points of each quadrilateral, either a 4-CST (Constant 

Strain Triangles) or an averaged 2-CST quadrilateral 

element, with respect to the global system (Figure 38) 

AX ( n ) , BX( n ) , CX(n>, DX(n> 

AY(n ) , BY(n ) , CY(n ) , DY(n ) 

where n * 1, 2, 3. N, and N Is the number of the total 

elements. Thus, we require Ircore storage with 8 x N words 

for the eight variables. Actually, we need only assign the 
relative coordinates of one elment since we are temporarily 
dealing with regular sized elements. The Introduction of 
the relative coordinates Is based on the fact that both the 
shape function END and the strain-displacement matrix [B1 
of the trlanglular elements depend only on the relative 

coordinates of the element. Thus, In the new version of 

AFEA, we need only eight Incore words for these eight 

variables. 
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(4) RESTART/BACK-UP OPTION 

In the new version of the two-dimensional computer 
code, we have Included a restart/back-up option. That Is, 
we save the displacements and velocities of every nodal 
point after a certain number of time steps. The new back-up 
results will overwrite on the same tape to replace the old 
back-up results. If the Integer IFLAG has the value zero, 
the option will start from the beginning; any other values 
of IFLAG will cause the option to restart from Its last 
executed time step. This feature greatly facilitates 
efficiency of computation In case of a computer failure. 




(5) AVERAGING 2CST FORMULATION 


In the old version of our two-dimensional finite 
element codes, <SH-wave and Elastic Wave Cases>, we used 
4-CST (Constant Strain Triangles) quadr11atera1 elements In 
order to avoid the space-grid skewness. The degrees of 
freedom of the Internal node, such as node 5 In Figure 39-A, 
are statically condensed. The degrees of freedom of such 
nodes do not appear In the global assemblage equations so 
that the storage during the computation Is reduced 
considerably. The 4-CST formulation has been a very good 
approach to obtain solutions for static cases. However, in 
the 4-CST algorithm for the dynamic case, It was assumed 
that the Internal condensed nodal points for each 
quadrilateral element were subjected to no Inertial forces. 
In the present modified version of these codes, an averaged 
2-CST algorithm Is adopted. In Figure 39-B, a quadrilateral 
is composed of two triangles. Two distinct ways of 
subdividing the quadr11atera1 are given. Therefore, the 
assembled quadrilateral stiffness matrix can be obtained by 
a linear Interpolation without any Internal condensation. 
Based on the averaged 2-CST formulation, a better- 
approximation of the solution can be obtained. 
















(6) Noda 1 -Po 1 nt-Or 1 ented Approach 

In the old version of our two-dimensional finite 
element codes, the global stiffness matrix Is obta-fned by 
the assemblage of the stiffness matrix of each Individual 
quadrilateral element for each time step In the time 
1ntegrat1 on. 

In the present modified version of the codes, the 
global matrix Involves only the non-zero stiffness matrix 
members for each nodal point before time Integration. 

The advantages of the old codes are: 

(1) The Incore storage required Is very little. Only 
the products of the global stiffness matrix and 
displacements are stored. 

(11) There will be less chance of making mistakes, 
particularly for source mechanism problems. 

The advantages of the new codes are: 

(1) The computing time Is drastically reduced since 
local assemblages for each time step are avoided. 
The amount of the reduced time depends upon the 
size, and the property of each individual problem. 
For a model of 40x40 elements, the computing time 
using the new version program Is only 1/5 of that 
using the old version. 

(11) The disk storage for the problem with Irregular 
element meshes Is the same as that of the 
problem with a regulai size element mesh. 

In both the versions, the 2-CST formulation and 
effective excitation method are used. The restart/back-up 
option Is also Introduced Into the new version codes. 











IV. FINITE ELEMENT SOURCE MECHANISM 


The finite element simulation of a non-directional line 
or a point source, ( cylindrical source for two-dimensional 
problems, or a spherical source for three-dimensional 
problems), for the two- or three-degree-of-freedom elastic 
wave case has been a very challenging subject. The simplest 
or the most natural way to excite the elastic medium is to 
apply a directional forcing function to the structure. 
Likewise, for the one-degree-of-freedom problem, special 
effort Is needed to simulate a coupled source. By using an 
energy-shar1ng-noda1-po1nts technique, we have successfully 
simulated different types of sources. We have studied the 
two-dimensional source mechanism problems for the cases of 
SH-waves and elastic waves (both P and S), respectively. 
For the SH-wave case, two type of sources are considered: 
(a) a concentrated line source, and (b) a concentrated 
coupled lino source. For the elastic wave case, three types 
of source are considered, viz., (a) a d1rect1ona1-force line 
source, (b) an omni-directional line source, (c) a 
s1ng1 e-coup 1e-w1thout-moment line source, (d) a 
s1ng1 e-coup 1e-wlth-moment line source, and (e) a 
douple-couo1e line source. 

Consider a two-dimensional whole space elastic medium, 
with Vp * 5,800 m/sec, v s * 2,900 m/sec, ^ *2.67 gm/cm^ . 
Figure 40 shows the 30x30 finite element mesh. Point 
S10,0) is the location of the applied source. Using 


this model, we have studied the following cases: 
































(1) ONE DEGREE OF FREEDOM SH-WAVE FROBLEM 

(a) Concerrtrated Line Source 

As we mentioned previously, this is the mosrt natural 
way to apply an external source. At Point S(0,0) In Figure 
40, a forcing function, F^(t) as given In Figure 2, is 
applied. Figures 41-A to 41-H are the displacement fields 
at 0.1 sec Intervals. The cy11ndr i ca1ly spreading wave can 


be observed. 

















(11) Concentrated Coupled Line Source 

In order to simulate a coupled line source at point 
S(0,0) of Figure 40, we let the points S^(0,+0), and 
S2<0,-0) occupy the same location as point S and apply the 
forcing function F^(t) to each Individual nodal point but 
In the opposite directions as shown In Figure 42. In other 
words, the displacement fields are w-^F.+fl.t) and 
w^(0,-0,t) at S and , respectively. At the source 
points In Figure 42, dot means outward from the plane of the 
paper, cross means Inward Into the plane of the paper. 
Points A, B, C, and D are the nodal points surrounding point 
S (Figure 40). Figure 43-A to 43-H are the snapshots of the 
synthetic displacement field w(x,y,t) at 0.1 sec Intervals, 
showing the radiation pattern of a coupled source. 











































microcopy resolution test chart 
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(2) TWO DEGRE£ OF FREEDOM ELASTIC WAVE PROBLEM 


(a) S1ng1e-D1rectIona 1 Force Line Source 

For the elastic wave case,, one of the most: common 
sources Is the directional line source. In Figure 44, a 
vertical line source, F Q <t) as given In Figure 18, is 
located at point S. Figures 45-A to 45-H are the snapshots 
of the displacement components u and u . Figures 46-A to 

x y 

46-1 are the snapshots of the absolute values of the 
displacement resultants. The time Interval between 
snapshots Is 0.0575 sec. 
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the Detonation of the Directional line Source 
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(b) Omni-Directional Force Line Source 


In order to simulate an omni-directional line source at 
point S(0,0) In Figure 40, we let the points S^(-0,+0), 
S2< + 0, + 0>, S 3 <+0,-0), and S,j(-0,-0) reside at the same 
location as point S and apply a forcing function, F 0 (t), to 
each Individual point In the manner shown In Figure 47. 
Figures 48-A to 48-H are the snapshots of the displacement 

components u and u . The time Interval between snapshots 

x y 

Is also 0.0575 sec. The cy11ndr1ca1ly spreading radiation 
patterns of the displacement magnitude can be observed in 
Figure 49-A to 49-0. 
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(c) Couple-WIthout-Moment Line Source 


Figure 5 S3 shows -the mechanism of a 
couple-wlthout-moment line source -for the elastic wave case. 
Points S^(+0,0) and S 2 <occupy the same location and 
are subjected to two forcing functions, F (t>, having the 
same magnitude but In the opposite directions, creating a 
couple without moment. The radiation patterns of the 

displacement components u and u are given In Figures 

x y 

51-A to 51-B. 










































(d> Couple-W1th—Moment Line Source 

Figure 52 shows tire mechanism of a 
couple-wlth-moment line source for the elastic wave case. 
Points and occupy the same location and 
are subjected to a moment forcing function F Q (t), with the 
same magnitude but In the opposite directions. The 
radiation patterns of the displacement components u^ and 
u are given In Figures 53-A to 53-B. 




















J 




.. *vf 



« 

■ mm 





«WWMWHMlk.S; . . *, I 

i^illffillHlIMlM .t I 

L •*- ;■ , . 

4 ^^: 7 iSglj | 



% r> 


















(d) uoub1e-Couple Force Line Source 


In order to simulate a double-couple force line source 
at point S(0,0) In Figure 40, we let the points Si(-0,*0), 
S 2 (+0,+0), S 3 < ■*■0,-0), and S^I-0,-0) at the same location as 
point S and apply a forcing function F Q (t> to each 
Individual point In a manner shown In Figure 54. Figures 
55-A to 55-H are the snapshots of the displacement 
components u x and Uy. The time Interval between snapshots 
Is also 0.0575 sec. 














































< 3) REMARKS ON THE FINITE ELEMENT SOLUTIONS 


(a) In the preceding seven cases, the non-reflecting 
boundaries techniques are used on all the four boundaries. 
Without using the non-ref1ec11ng boundaries, undesired 
reflections would appear at the 15th time step after the 
detonation of the source, or 0.3 seconds for the SH-wave 
cases and 0.1725 seconds for the models of the elastic wave 
cases . 


(b) All the snapshots of the displacement field for the 
preceding five elastic solid cases are presented In a 
rectangular coordinate system. It may provide better 
physical understanding If we transform the results into the 
polar system (r,8). Figures 56 show the radiation patterns 
for the radial displacement u r and the transverse 
displacement u t due to the five different types of source 
(Figure 56-A). The snapshots shown In Figure 56-B are the 
displacement fields at the 30th time step, 0.345 seconds 
after the detonation of the source. These results 
qualitatively agree with the classic analytical solutions of 
the three-dimensional source mechanism given by Honda 
( 1952 ) . 
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(B) Omni-Directional Force. 
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(E) Double-Couple Force. 


Figure 56 - A. Five Different Types of Source. 
















<c) Honda (1952) showed that the combination of two 
coplanar-couple forces with moment In the opposite 
rotational direction (Figures 57-A and 57-B) Is equivalent 
to that of two coplanar-couple forces without moment. (or a 
comlnatlon of compression and tension, which are 90 degrees 
out of phase (Figure 55-C). In the finite element results, 
Figures 58-A and 58-B are the snapshots of the radial and 
transverse components of the displacement field, u r and u t , 
associated with the cases of Figure 57-A and 57-B. Figure 
58-C Is the combined results of 58-A and 58-B. These 
combined results are exactly the same as the displacement 
field due to the case of Figure 55-C. The snapshots are 
taken at the 25th time step, 0.2875 seconds. 
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( 4 ) DEMONSTRATION OF THE WAVEFORMS STRONGLY DEPENDS ON THE 
NATURE OF SOURCE 

It is well known that the elastic response due to a 
spherical pressure pulse Is dominated by P waves, while the 
elastic response due to a directional point force Is rich In 
shear and Rayleigh waves In the near field. In order to 
demonstrate this difference In waveforms due to two 
different types of source, an example of an anticlinal model 
Is given In Appendix A (A part of the work In Appendix A Is 
supported by PROJECT MIDAS). 
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V. THE TWO-DIMENSIONAL FINITE ELEMENT COMPUTER PROGRAM FOR 
ELASTODYNAMIC PROBLEMS <THE USER’S GUIDE) 

The program described here Is specifically coded for solving 
a class of elastodynam1c transient problems in elastic wave 
scattering and diffraction for the two-dimensional case. This 
element-oriented finite element code, designated as ALFE2.A, has 
been adapted to the VAX 11/780 computer at AFGL last October 
( 1984 ). 

In the computer code ALFE2.A, we have addressed our 
particular attention toward the areas of: (a) practical 
difficulties of using enormously large in-core storage without 
Increasing I/O time, and (b) solving sparse matrix without 
calculating zero-matrix elements to save computing time. We 
succeeded In overcoming these serious difficulties to a certain 
degree. 

The present computer code ALFE2.A Is written In standard 
FORTRAN IV language and has been executed successfully on the 
PRIME 750 computer. 

(1 ) FORMULATION 

The general procedure for using the finite element method In 
solving e1astodynam1c problems In continuum mechanics Includes 
Idealization of the problem, finite element space discretization, 
formulation of the set of simultaneous equations to be solved, and 
time Integration. We summarize the ALFE2.A code as follows: 






(a) Basic Formulation 


Based on the principle of virtual work, the equations of 
motion for e1 astodynam1cs are: 

0] \“] ♦ l>]f^ - 

where juj, juj are the displacement and acceleration vectors of the 
finite element assemblage, and [m] Is the assembled mass matrix, 
[k] the assembled stiffness matrix, and |V(t)j Is the 
applied time-dependent load vector. 

(b> Spatial Elements 

For the present version of Program ALFE2.A, the 
two-d(mensfona1 solid unit element In the nodal numbering system 
Is that of an four-node solid quadrilateral element. Each 
quadrilateral element Is comprised of two constant strain 
t r 1 a r.g 1 < 2CST ) . Two distinct ways of subdivision are used In a 

manner as shown In Figure 59. A standard right-handed Cartesian 
coordinate system Is used for the global coordinates. 

(c) Interpolation Functions for a Triangle 

The functions of the displacement field 
approximated by two linear polynomials 
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The *1x constants ^ (1 ■ 1,2,3) can be evaluated by 
solving two sets of three simultaneous linear equations. These 
simultaneous equations can be obtained by equating the values of 
the displacements at the triangle nodes. 1, J, and k (Figure 6JBT). 
The coefficients . and ^ can be solved In terms of the nodal 
displacements and nodal coordinates. Once these coefficients are 
determined, we can write the displacement of an arbitrary point 
within the triangle element as 



where I Is a 2x2 unit matrix, N., N., N^, or ['n') are the shape 
functions, which are expressed In terms of the nodal 
coordinates, ^u^ Is a column vector with six nodal displacement 
components of the triangle element. The strain-displacement 
matrix [b] can bo defined explicitly In terms of ["nJ. 

(d) Triangle Stiffness Matrix 

Construct the elastic matrix CD] and carry out the matrix 
multiplications CB] T CD3CB]. After Integration of CB1 T CD][B] over 
the area, the triangle area. In this linear case, a 6x6 triangle 
stiffness matrix is obtained. 

The Young's Modulus and Poisson's ratio (or the Lame's 
constants) of the material are required to define the stiffness 
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matrix of the solid element. 

(e) Quadrangle Stiffness Matrix (Averaged 2CST> 

The final 8x8 element stl'ffness matrix Is obtained by 
assembling two 6x6 stiffness matrices. No condensations of 
Internal degrees of freedom are required, since there are no 
Internal nodes In this assemblage. 

(f) Mass Matrix 

By assuming that all the mass attributed to a node Is 
concentrated at that node, the mass of the quadrilateral 
element Is lumped welghtedly to the four nodal points. The 
mass density of the material Is required to define the mass 
matrix of the element. 

(g) Global Stiffness and Mass Matrices Generation 

Assemble the stiffness and mass matrices of the 
Individual quadrilateral elements with a common nodal point 
to obtain the equations of motion at the nodal point. A 
final banded characteristic global stiffness matrix of order 
Ixl Is obtained, where I ■ NxNF, N Is the number of total 
nodal points, and NF Is the degrees of freedom of each nodal 
point. The banded global stiffness matrix Is highly sparse. 
The mass matrix Is diagonal. 
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(h) External Load Vector and Prescribed Boundary Conditions 

The external forcing function F(t) may be read in as a 
single vector. Traction boundary conditions are Incorporated Into 
the external load vector automat leally,1.e. for a free boundary, 
the load terms associated with the nodes on the boundary are zero. 
For specified displacements, the stiffness matrix has to be 
mod 1f1ed. 

Any nodal point may be chosen to be one of the four types of 
boundary conditions, depending on the code number of each nodal 
point. 

NZRPT12*N-1 )■ 1 -- the nodal point Is free to traction in 

the x-dlrectlon 

NZRPT(2*N) «1 -- the nodal point Is free to traction In 

the y-dlrectlon 

NZRPT(2*N-1)-0 -- the nodal point Is fixed In the x-dlrectlon 
NZRPT(2*N) -- the nodal point Is fixed In the y-dlrectlon 

The character N Is the array number of the nodal point. 

(1) Time Integration Scheme 

The following time integration scheme 

|u(t+At)j= ju(t)j ♦ |u(t)^At 

ju ( t+At) j = ju(t)j + [M] j ju ( t+4t) + [m] 1 1 1 )| A t 

Is carried out step by step to determine the displacements and 
velocities In the x, and y directions. The quantity £t Is the 
time step Increment. 
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(2) PROGRAM DI SCR IPTION 

The ALFE2.A code consists of Main Program and seven 
subroutines. These subroutines carry out most of computational 
steps for spatial stiffness and mass matrix assemblage, and are 
controlled by Main Program. Following are descriptions of the 
functions of Main Program and subroutines. 

MAIN P7CC7.AM: Solves overall time-history response of the ent 1 re 

finite element model, and prints the displacements, 
velocities and accelerations. 

SUBROUTINES s 

DATAIN - reads In the Input data. An automatic finite element 
mesh generator is Introduced. 

TRIANG - formulates the stiffness matrix of a single triangle 
e1ement. 

QUAD - formulates the constitutive relations. 

STIFT - formulates the stiffness and mass matrices of a 
single quadr11atera1 element assembled by two 
triangular elements. 

FRMSTF - formulates the global stiffness and mass matrices, 
finds the products of the global stiffness matrix 
and displacements. 








FORCE - Introduces the forcing function. 

NEW - formulates the quadr11atera1 stiffness matrix by 
an alternative approach. 

The flow charts of Main Program ALFE2.A and all the 
subroutines are appended at the end of this section. 

Main Program first calls subroutine DATAIN to read and echo 
print the Input master control data, and material properties for 
each element. 

The coordinates for each nodal point and tht numbering of a 
system of nodal points Interconnection are generated 
automatically. As the elemental data are Identified, subroutines 
FRMSTF, STIFT, QUAD, TRIANG are then called to form the global 
stiffness matrix and global mass vector matrix based on the 
quadrilateral stiffness and mass matrices. The products of the 
global stiffness matrix and displacement are computed at the same 
time as the global stiffness matrix Is assembled. 

Subroutine FORCE Is called to formulate the applied load 
vector matrix. Time Integration Is then carried out to determine 
the nodal displacement and velocity components In each time step 
1nc rement. 
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(3) MULTILAYETTED MESH GENERATION 

The computer Code ALFE2.A Itself Is applicable to any 
two-dimensional geological and stratlgraph1c structures- In 
this preliminary version ALFE2.A, the regular-size finite 
element meshes are generated by an automatic mesh generator. 
The basic element or zone of the generalized grid Is a 
rectangle. Figure 61 shows a 1-2x 12 diagram. 

Figure 62 shows the grid-point numbering convention for 
a 3x3 mesh, l.e. NROW * {MCOL-1 ) ■ 3, where 

(MCOL-1) ■ number of elements In the x-dlrectlon. 


NROW 


number of elements In the y-dlrectlon. 


The grid points are numbered In the order: Top to Bottom , 
Right to Left . 

Figure 63 provides an example showing the convention 
for numbering the elements for 3x3 mesh. The elements are 
numbered In the order Top to Bottom , and Right to Left . In 
this mesh generator, (MCOL-1), and NROW are not necessarily 
to be the same. 


In this version of ALFE2.A, the automatic finite 


element mesh generation subroutine Is limited to nine 
horizontal layers. That Is, the finite element mesh can 
model various different combinations corresponding to one to 
nine different materials. For example. In lines 3-11 of the 
data Input, If we assign the same elastic parameter 
constants for each lijyer, the mesh will simulate a 
homogeneous medium. If we assign different elastic parameter 
constants to the first layer from those of the other eight 
layers , then we obtain a one-layered medium, and so forth. 































F igure 63. An Example Shoving the Convention for 
Number inq of Elements in PROGRAM ALFE2.A. 













(4) GUIDE FOR DATA INPUT 


BASIC 

PARAMETERS (FORMAT 

315 ) 


COLUMN 

VARIABLE 


DEFINITION 

1-5 

I FLAG 


Back-up Option 

IFLAG-0 : Start from beginning 
IFLAGP0 s Restart from Its last 
executed time step 

6-10 

NBSTEP 


Back-up by NBSTEP time steps 

11-15 

I SRC 


Code for source type,(temporar1 
set to "l", means directional 
force source) 

BASIC 

PARAMETERS (FORMAT 

215, 

F15.8 > 

COLUMN 

VARIABLE 


DEFINITION 

1-5 

NMATRL 


Number of materials 

6-10 

NSTEP 


Total number of time steps 
to be calculated 

11-25 

DELTAT 


Time-step Increment 


MATERIAL PROPERTIES (FORMAT 3F15.4) 


COLUMN 

VARIABLE 

DEFINITION 

1-15 

VP 

Compress 1onaT wave velocity 
of the material 

16-30 

VS 

Shear wave velocity of the 
mater 1 a 1 

31-45 

RHO 

Mass density of the material 
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LOADING POINT 

AND MATERIAL 

-ELEMENT CONNECTION (FORMAT 1115) 

COLUMN 

VARIABLE 

DEFINITION 

1-5 

NS 

Force applied on NS-th 
vertical nodal point 

6-10 

NROWSO 

Force applied on NRORSO-th 
horizontal nodal point 

11-15 

MH1 

Vertical number of elements 
for the first layer 

16-20 

MH2 

Vertical number of elements 
for the first two layers 

21-25 

MH3 

Vertical number of elements 
for the first three layers 

26-30 

MH4 

Vertical number of elements 
for the first four layers 

31-35 

MH5 

Vertical number of elements 
for the first five layers 

36-40 

MH6 

Vertical number of elements 
for the first six layers 

41-45 

MH7 

Vertical number of elements 
for the first seven layers 

46-50 

MH8 

Vertical number of elements 
for the first eight layers 

51-55 

MH9 

Vertical number of elements 
for the first nine layers 

COORDINATE AND 

NODAL POINT 

CONTROL (FORMAT 2I5.2E16.4) 

COLUMN 

VARIABLE 

DEFINITION 

1-5 

NROW 

Number of element In the 
y-d1rec11 on 

6-10 

MCOL 

Number of nodal points In 
the x-dlrectlon 

11-15 

XX 

Size of quadrilateral elment 
In the x-dlrectlon 

16-20 

YY 

Size of quadrilateral element 
In the y-dlrectlon 
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(5) DEFINITION OF OTHER MAJOR VARIABLES 


VARIABLE DEFINITION 


NODEPT 

NELMNT 

NF 

MAXNOD 

NEQTNS 

NELMAR 

AX. BX.CX.DX, 

AY, BY.CY.DY 

S(6,6) 

AREA 

STIFMX(8,8) 
SMX(8,8) 

E< N ) 

PNU(N ) 

P L UMP(N ) 

PMU <N) 

D(3,3 ) 

B(3,6 ) 

NZRPT 

UUU 

VVV 

UUI 

VVI 

XXX 

MBAND 


Number of the totel nodal points 

Number of the total elements 

Degrees of freedom of each nodal point 

Maximum number of nodal points for each 
element (Quadrangle) 

Total number of equations 

Element connecting array 

Coordinates of the four nodal points of 
quadrilateral element with respect to 
the global system 

Stiffness matrix of triangular element 

Area of triangular element 

Stiffness matrix of quadrilateral element 

Mass matrix of quadrilateral element 

Young's Modulus of the N-th material 

Poisson's ratio of the N-th material 

Lame's constant > v of the N-th material 

Lame's constant yU of the N-th material 

Constitutive matrix CD] 

Stra1n-d1sp1acement matrix CB3 

Boundary condition code 

Displacement component of each nodal 
po 1 nt 

Velocity component of each nodal 
point 

Displacement component of each nodal 
point for different type of non¬ 
reflecting boundaries 

Velocity component of each nodal 
point for different type of non¬ 
reflecting boundaries 

Matrix product of CK] juj 

Half-band width of assemblage equations 
(Irrelevant to the present code ALFE2.A) 








































SUBROUTINE DATAIN 






































SUBROUTINE STIFT 
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APPENDIX 

AN ANTICLINAL ELASTIC MODEL 

(Demonstration of the Radiation Patterns Strongly 
Depends on the Nature of Source) 

In order to demonstrate the differences of a vertical 
forcing line source and an omni-directional line source, we 
consider a two-dimensional elastic anticlinal model as shown 
In Figure 1-A with the compress 1ona1 velocities of each 
medium: v ^ *7 ,000 ft/sec, Vp 2 *5,900 ft/sec, and 

v p3 ”6,400 ft/sec, respectively, and the densities 

*2.1 gm/cm 3 , p 2 ”1*8 gm/cm 3 , and *1.9 gm/cm 3 . The 

Poisson's ratios are assumed to be 0.25 for each medium. The 
horizontal range Is 2,000 ft, the depth to the top of the 
anticline Is 1,000 ft. Two different kinds of source are 
applled to point S, namely, < 1 ) a vertical force, and (11) an 
omni-directional force. The source Is located at a depth of 
200 ft at the center of the horizontal range. The source 
function Is taken to be the first derivative of Gaussian with 
a duration of 0.04928 sec. The center frequency of this 
source function Is approximately 41 Hz. The magnitudes of all 
the forcing functions for a vertical and an omni-directional 


11 




































line source used are the same. Taking the advantage of 
symmetry with respect to the z-axls, we use only one half of 
the model (Figure 1 — B >. The finite element mesh Is 
automatically generated, and Is given In Figure 1-C. The 
smallest element size , As ,1s 10 ft, and the time step ,At, 
Is taken to be 0.88 msec. Three synthetic seismograms for 
each case, one along the free surface, and two along the 
vertical planes OA and CD (see Figure i-B)are recorded. 
Figures 2-A and 2-B are the vertical and hor1zonta1-component 
displacements along the free surface for the vertical force 
line source problem. There are 101 traces plotted. Figures 
3-A,B and 4-A,B are the vertical and horizontal components of 
the displacement along the vertical planes OA and CD for the 
vertical-force source problem, respectively. In Figures 3 and 
4, there are 71 uneven spaced traces plotted, according to the 
actual finite element mesh locations. Figure 5, 6, and 7 are 
the vertical and hor1zonta1-components of the displacements 
along the free surface, the vertical planes OA thed CD, 
respectively for an omni-directional line source problem. 

For 1nterpretat1 on, we adopt the following notations for 
the arrivals In Figures 2 to 7: 

(1) subscripts denote medium, prime denotes the last 


med1 urn. 
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(2) the arrival of compreysl orra.1 wave denotes P and that 
of shear wave S. 

(3) the wave arrival from the source to the free surface 
Is In 1.c. . 

(4) Reflected P once from source to surface, Is denoted 
by P^P^. reflected again from the Interface of-med1 urn I 
and II Is p^P^P^ and 80 on • 

In Figures 2 and 5, the primary arrivals, p^ , 
s^P^-dIffracted, S, Rayleigh, and the mu 111 - ref 1ected P waves 
from the bottom of medium I, labeled as P i P i* P 1 P 1 P 1 P 1’ 
P 1 P 1 P 1 P 1 P 1 P 1’ have been identified. It Is clearly seen that 
the responses of the free surface, and along the vertical 
planes OA and CD are significantly different between a 
vertical force line source and an omni-directional force line 
source. In Figure 2, the vertical displacement Is domlnanted 
by Rayleigh waves, which also are major contrlbuters to the 
horizontal displacement. The s -^P -^-d Iff racted waves are masked 
by the arrivals of multiple P reflections. The 
s^P^-dIffracted waves start on the free surface around 150 ft 
away from point 0 with an approximate arrival time of 0.0618 
sec after the Initial detonation. In Figure 5, there are 
virtually no S waves arrivals before the waves reach any 
boundary. The predominant wave In this case Is P wave. 
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Comparison of Figure 3-A with Figure 6-A, the 
vert 1ca1-components of the displacement along the plane OA due 
to a vert 1ca1-force line source (VFS) and an omn1-d1rect1ona1 
line source (OLS), respectively shows that: 

<1> not only wave forms are different, but also p^ and P^ 
are stronger for VFS than these for OLS. 

(2) the late arrivals of F^ and p^P^P^R^ . are 

clearer for OLS than these for VFS. 

(3) except the wave forms are different, aro 

clearly 1 dent 1f11ab1e for both VFS and OLS. 

I 

(4) the Initial displacements of p^^P^P-jP^ are stron 9 er 
for VFS than these for OLS 

(5) since the plane OA contains the source S, there are 
no displacements In the hor1zonta1-components. 

Comparison of Figure 4-A with Figure 7-A, the 
vert 1ca 1-components of the displacement along the vertical 
plane CD due to VSF, and OLS, respectively shows that: 

(1) except shifted In the arrival times and smaller in 
magnitudes, the characteristics of the primary arrivals namely 
p^ and Pare virtually Identical to these along the plane 
OA. 
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Figure 5 -B. synthetic seismogram along free surface 



































































































































































FIGURE a - B . SYNTHETIC SEIIHOORRH ALONG VERTICAL PLANE OA 
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(2) because of CD Is very close to OA, the vertical 
displacements of S waves are small and not clearly 

Identifiable In either Figures 4-A or 7-A. 

Comparison of Figure 4-B with Figure 7-B, the 

horizontal-components of the displacement due to VFS and OLS 
shows that: 

(1) S waves are easily Identifiable for both VFS and OLS. 

(2) S waves are stronger for VFS than these for OLS. 

(3) p^S-^ and p^P-^ are Identifiable for OLS and only p^P^ 
Is identifiable for VFS. 

(4) Mutlples, such as, p i I> i P 2 P 3 P 1 . are qu1te 

clearly recorded for both VFS and OLS. 


CONCLUSION 


From the above results of the simulation of VFS and OLS, 
It Is evident that seismic radiation patterns are strongly 


dependent on the nature of the source. 
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